version 8.0 
log using logfilename.log, replace
#delimit ;

*     ***************************************************************** *;
*     ***************************************************************** *;
*       File-Name:      limited.do                                      *;
*       Date:           Month/Day/Year                                  *;
*       Author:         MRG                                             *;
*       Purpose:        Do-file for creating figures showing the        *;
*                       marginal effect of X on Y with a single         *;
*                       modifying variable where the dependent variable *;
*                       is limited.  The estimated model is:            *;
*                       Y = b0 + b1X + b2Z + b3XZ + epsilon             *;
*       Input File:     name_of_datafile.dta                            *;
*       Output File:    logfilename.log                                 *;
*       Data Output:    None                                            *;
*       Previous file:  None                                            *;
*       Machine:        Office                                          *;
*     ****************************************************************  *;
*     ****************************************************************  *;

use h:\name_of_datafile.dta;


*     ****************************************************************  *;
*                           Estimate Probit Model                       *;
*      Y* = b0 + b1X + b2Z + b3XZ + b4Control1 + b5Control2 + epsilon   *;
*     ****************************************************************  *;

probit Y   X   Z   XZ  Control1 Control2;

*     ****************************************************************  *;
*       Take 10,000 draws from the estimated coefficient vector and     *;
*       variance-covariance matrix.                                     *;
*     ****************************************************************  *;

preserve;
drawnorm MG_b1-MG_b6, n(10000) means(e(b)) cov(e(V)) clear;

*     ****************************************************************  *;
*       To calculate the desired quantities of interest we need to set  *;
*       up a loop.  This is what we do here.                            *;
*     ****************************************************************  *;
*       First, specify what you quantities should be saved and what     *;
*       these quantities should be called.                              *;
*     ****************************************************************  *;

postutil clear;
postfile mypost prob_hat0 lo0 hi0 prob_hat1 lo1 hi1 diff_hat diff_lo diff_hi 
            using sim , replace;
            noisily display "start";
            
*     ****************************************************************  *;
*       Start loop.  Let `a' be the modifying variable Z and let this   *;
*       run from min to max in the desired increments.                  *;
*     ****************************************************************  *;                                  

local a=0 ;
while `a' <= 35 { ;

    {;

scalar h_X=5;
scalar h_Control1=2;
scalar h_Control2=3.4;
scalar h_constant=1;

    generate x_betahat0 = MG_b1*h_X
                            + MG_b2*(`a')
                            + MG_b3*h_X*(`a')
                            + MG_b4*h_Control1
                            + MG_b5*h_Control2
                            + MG_b6*h_constant;
    
    
    generate x_betahat1 = MG_b1*(h_X+1)
                            + MG_b2*`a'
                            + MG_b3*(h_X+1)*(`a')
                            + MG_b4*h_Control1
                            + MG_b5*h_Control2
                            + MG_b6*h_constant; 
    
    gen prob0=norm(x_betahat0);
    gen prob1=norm(x_betahat1);
    gen diff=prob1-prob0;
    
    egen probhat0=mean(prob0);
    egen probhat1=mean(prob1);
    egen diffhat=mean(diff);
    
    tempname prob_hat0 lo0 hi0 prob_hat1 lo1 hi1 diff_hat diff_lo diff_hi ;  

    _pctile prob0, p(2.5,97.5) ;
    scalar `lo0' = r(r1);
    scalar `hi0' = r(r2);  
    
    _pctile prob1, p(2.5,97.5);
    scalar `lo1'= r(r1);
    scalar `hi1'= r(r2);  
    
    _pctile diff, p(2.5,97.5);
    scalar `diff_lo'= r(r1);
    scalar `diff_hi'= r(r2);  

   
    scalar `prob_hat0'=probhat0;
    scalar `prob_hat1'=probhat1;
    scalar `diff_hat'=diffhat;
    
    post mypost (`prob_hat0') (`lo0') (`hi0') (`prob_hat1') (`lo1') (`hi1') 
                (`diff_hat') (`diff_lo') (`diff_hi') ;
    };      
    drop x_betahat0 x_betahat1 prob0 prob1 diff probhat0 probhat1 diffhat ;
    local a=`a'+ 1 ;
    display "." _c;
    } ;

display "";

postclose mypost;

*     ****************************************************************  *;                                  
*       Call on posted quantities of interest                           *;
*     ****************************************************************  *;                                  

use sim, clear;

gen MV = _n-1;

graph twoway line diff_hat MV, clwidth(medium) clcolor(blue) clcolor(black)
        ||   line diff_lo  MV, clpattern(dash) clwidth(thin) clcolor(black)
        ||   line diff_hi  MV, clpattern(dash) clwidth(thin) clcolor(black)
        ||  ,   
            xlabel(0 10 20 30, labsize(3)) 
            ylabel(-.1 0 0.1 0.2, labsize(3))
            yscale(noline)
            xscale(noline)
            yline(0)
            legend(off)
            title("Marginal Effect of X on Y As Z Changes", size(4))
            subtitle(" " "Dependent Variable: Y" " ", size(3))
            xtitle(Z, size(3.5)  )
            ytitle("Marginal Effect of X", size(3.5))
            xsca(titlegap(2))
            ysca(titlegap(2))
            scheme(s2mono) graphregion(fcolor(white));

*     ****************************************************************  *;
*                 Figure can be saved in a variety of formats.          *;
*     ****************************************************************  *; 

graph export  h:\figure1.eps, replace;

translate @Graph h:\figure1.wmf; 

*     ****************************************************************  *;
*                                   THE END                             *;
*     ****************************************************************  *;


exit; 
